Ultrasound imaging has got immense importance as diagnostic tool in medical applications for its low cost and non-invasive imaging modality. But the resolutions and speckle noises are the dominant issues, which reduces its utility in some medical diagnostics applications. Many speckle reduction techniques has already been proposed so far. All these techniques are applied either on the raw scan data in pre-processing stage (before scan conversion) or on the scan-converted image as post-processing operation. It is found that the image quality is relatively better if filtering is applied in the preprocessing stage (before scan conversion) rather than post processing stage (after scan conversion). But in this case, the amount of data to be handled is larger. This is true for all the popular types of speckle reduction filters. Furthermore, after noise reduction from raw data, interpolation is performed as part of scan conversion.
Ultrasound imaging modality is one of the most widely used imaging modality in diagnostic medical applications because, it is noninvasive, non-ionizing, real-time, practically harmless to human body, portable and cost effective. Unfortunately, the image quality of medical ultrasound imaging system is limited by some physical phenomena underlying in the acquisition system.
Speckle noise generated in the image is one such limitation. Ultrasound speckle is a quasi-random phenomenon as discloses in Czerwinski, R. N., Jones, D. L., William D. O'Brien, Jr., “Ultrasound Speckle Reduction by Directional Median Filtering”, Proceedings, International Conference on Image Processing, Vol: 1, (1995) which occurs due to backscatters ultrasound pulses from the rough surface of the internal soft tissues. Thus ultrasound speckle is similar in origin to laser or radar speckle. It degrades the resolutions, contrast and obscures the underlying anatomy and makes human interpretation and computer-assisted detection techniques difficult and inconsistent as disclosed in Michailovich Oleg V. Tannenbaum Allen, “Despeckling of Medical Ultrasound Images”, IEEE Trans. on Ultrasonics Ferroelectrics and Frequency Control, Vol. 53, No. 1, pp. January (2006).
Hence, reduction of speckles is one of the most important challenges to the ultrasound system designers'. Many attempts are made by the engineers and scientists to develop speckle reduction methods during last three decades, and, many techniques have also been developed as taught in Vera Behar, Dan Adam, Zvi Friedman, “A new method of spatial compounding imaging”, Ultrasonics 41, pp. 377-384, (2003), Pai-Chi Li and Mei-Ju Chen, “Strain Compounding: A New Approach for speckle reduction”, IEEE Trans on Ultrasonics Ferroelectrics and Frequency Control, Vol. 49, No. 1, January (2002), Jong-Sen Lee, “Digital Image Enhancement and Noise Filtering by Use of Local Statistics”, IEEE Trans. on Pattern Analysis And Machine Intelligence, Vol. 2, No. 2, March (1980), Gupta N, Swamy M. N. S., Plotkin E., “Despeckling of Medical Ultrasound Images Using Data and Rate Adaptive Lossy Compression”, IEEE Trans. on Medical Imaging, Vol. 24, No. 6, pp. 743-754, June (2005). These methods are basically applied either on the raw scan data in the preprocessing stage (i.e. before scan conversion) or on the scan-converted image in the post-processing stage (i.e. after scan conversion).
Basic theory Ultrasound Speckle and speckle statistics: The ultrasound B-scan imaging process is a result of a set of complicated physical phenomena such as absorption, reflection and coherent scattering of ultrasound pulse-echo signal from scattering medium. The back-scattered echo is received and used to display as image. The images, formed by such a process, involve granular structure called speckle. Basically, ultrasound speckle is generated from phasors' summation of coherent scatterings within the resolution cell as it is scanned through the phantom. This phenomenon can be treated geometrically as random walk of component phasors as disclosed in Robert F. Wagner, Stephen W. Smith, John M. Sandrik, H. Lopez, “Statistics of Speckle in Ultrasound B-Scans”, IEEE Trans. on Sonics and Ultrasonics, Vol. 30, No. 3, pp. 156-163, May (1983). If the number scatters within resolution cell is large, and the phase of the scattered waves is uniformly distributed within 0 and 2π independent of amplitude, the envelope of the complex phasor resulting from the summation of the scattered waves exhibits Rayleigh distribution.
The accumulation of the random scatterings can be represented by phasor summation of the scatterings as,
                    A        =                              ∑            i                    ⁢                                    a              i                        ⁢                          e                              j                ⁢                                                                  ⁢                                  φ                  i                                                                                        (        1        )            where each scatterer bears αi amount of signal and has a phase shift of φi. If αi and φi are assumed to be independent and identically distributed, the joint pdf of the real and imaginary component of the complex phasor can be given by central limit theory as,
                                          p                                          A                R                            ⁢                              A                I                                              ⁡                      (                                          A                R                            ,                              A                I                                      )                          =                              1                          2              ⁢              π              ⁢                                                          ⁢                              σ                2                                              ⁢                      e                          -                                                                    A                    R                    2                                    +                                      A                    I                    2                                                                    2                  ⁢                                      σ                    2                                                                                                          (        2        )            where σ2=E[AR2]=E[A12] is the second moments of the real and complex components AR and A1.
The envelope of the complex phasor can be calculated as,a=√{square root over (AR2+A12)}  (3)
Therefore the probability density function of the envelope is given by,
                                                                                          p                  ⁡                                      (                    a                    )                                                  =                                ⁢                                                      a                                          σ                      2                                                        ⁢                                      e                                          -                                                                        a                          2                                                                          2                          ⁢                                                      σ                            2                                                                                                                                                          ,                                                                        ⁢                              a                ≥                0                                                                                                        =                                ⁢                0                            ,                                                                        ⁢              otherwise                                                          (        4        )            
The function in equation (4) is known as Rayleigh pdf. The speckle pattern formed in the image under Rayleigh distribution is called “fully developed” pattern as disclosed in Dutt V. “Statistical Analysis of ultrasound Echo Envelope”, Ph.D. Thesis. Many other speckle statistics such as k-distribution, Rician distribution, Generalized gamma distribution, Weibull distribution, Nakagami distribution etc. are also considered in different literatures.
Most of the literatures of speckle reduction consider the multiplicative noise model for speckle noise as disclosed in Jain A. K., “Fundamentals of digital image processing” Book, Prentice-Hall, Inc and Kuan D. T., Sawchuk Alexander A. et al., “Adaptive restoration of images with speckle,” IEEE Trans. Acoustics, Speech and Sig. Proc., Vol. 35, pp. 373-383, March (1987). This multiplicative noise models for speckle is only a rough approximation, and ignore the correlation of the speckle that should be considered in speckle reduction.
Some important and popularly used speckle reduction techniques and the state-of-art briefly are mentioned below:
Speckle reduction techniques can be broadly categorized into three categories:                Compounding        Single scale spatial linear and nonlinear filtering        Multiscale method.        
Compounding techniques include                Spatial compounding        Frequency compounding        Strain compounding        
Underlying philosophy of compounding is the averaging of multiple images of the same target taken either from different positions, or with different frequencies or under different strains.
A number of works has been done on spatial compounding as taught in Fleming J. E. E., Hall A. J., “Two dimensional compound scanning-effects of maladjustment and calibration”, pp. 160-166, Ultrasonics, July (1968), Berson M., Roncin A., Pourcelot L., “Compound Scanning with an Electrically Steered Beam”, Ultrasonic Imaging 3, pp. 303-308, (1981) and Ping He, Kefu Xuet, Yiwei Wangt, “Effects of Spatial Compounding Upon Image Resolution”, Proceedings, 19th International Conference, IEEE/EMBS Oct. 30-Nov. 2, (1997) Chicago, Ill. USA.
In spatial compounding, multiple ultrasound images of a target are acquired by different spatial locations. Speckle in the common region of these images are partially correlated or not correlated. It is known that averaging of multiple images containing partially correlated or uncorrelated noises can reduce the effect of the noise. Hence, speckle can be reduced by forming a composite image averaging the acquired multiple images.
In the frequency compounding, the bandwidth of a radio-frequency (RF) signal is divided into a number of frequency sub-bands. Ultrasound signals from those bands are transmitted to form different images called sub-band images of the same target. A compounding image is, then produced by averaging the sub-band images. Speckles of the sub-band images are less correlated if the bandwidths of the sub-bands are narrower, since it is mainly determined by the difference of center frequencies, which is normalized by a 16 dB pulse envelope length of the sub-band signals as disclosed in Jin Ho Chang, Hyung Ham Kim, Jungwoo Lee, K. Kirk Shung, “Frequency compounded imaging with a high-frequency dual element transducer”, Ultrasonics 50, pp. 453-457, (2010).
Strain compounding as disclosed in Pai-Chi Li and Mei-Ju Chen, “Strain Compounding: A New Approach for speckle reduction”, IEEE Trans on Ultrasonics Ferroelectrics and Frequency Control, Vol. 49, No. 1, January (2002) exploits the decorrelation between signals under different strain states. Different strain states can be created using externally applied forces as the one used in sonoelastography. Such force produces three dimensional tissue motion. By correcting only in-plane motion, the images acquired under different strain states have similar characteristics except for speckle appearance caused by the uncorrelated out-of-plane motion. These images are combined for speckle reduction with less degradation in the in-plane spatial resolution.
However, all these compounding techniques suffers from different limitations:                They suffer from loss of temporal and/or spatial resolution.        Clinical application of strain compounding is limited        Contrast resolution for small objects may be degraded.        Complexity of the system increases.        
Single scale spatial linear and nonlinear filtering: Speckle reduction spatial filters perform smoothing according to local variance and local mean as discussed in different literatures i.e. Jong-Sen Lee, “Digital Image Enhancement and Noise Filtering by Use of Local Statistics”, IEEE Trans. on Pattern Analysis And Machine Intelligence, Vol. 2, No. 2, March (1980), Kuan D. T., Sawchuk Alexander A. et al., “Adaptive restoration of images with speckle,” IEEE Trans. Acoustics, Speech and Sig. Proc., Vol. 35, pp. 373-383, March (1987), Jong-Sen Lee, “Refined Filtering of Image Noise Using Local Statistics”, Computer Graphics And Image Processing 15, pp. 380-389, (1981), Jong-Sen Lee, “Speckle Analysis and Smoothing of Synthetic Aperture Radar Images”, Computer Graphics And Image Processing 17, pp. 24-32, (1981), Frost V. S., Stiles J. A., Shanmugan K. S., Holtzman J. C., “A model for radar images and its application to adaptive digital filtering for multiplicative noise,” IEEE Trans. on Pattern Analysis And Machine Intelligence, Vol-4, pp. 157-166, March (1982) and Bamber J. C., Daft C., “Adaptive filtering for reduction of speckle in ultrasonic pulse-echo images”, Ultrasonic, January (1986).
In the above mentioned filtering, techniques smoothing is increased in homogeneous region of the image and reduced or avoided elsewhere to preserve edges. These filters are basically adaptive filters. Adaptive filtering for reduction of speckle from ultrasonic pulse-echo images was proposed by Bamber J. C., Daft C., “Adaptive filtering for reduction of speckle in ultrasonic pulse-echo images”, Ultrasonic, January (1986). It proposed an adaptive two-dimensional filter which uses local features of image texture to recognize and maximally low-pass filter those parts of the image which correspond to fully developed speckle, while substantially preserving information associated with resolved-object structure. The filter is un-sharp masking filter and its output is mathematically given as,{circumflex over (x)}=x+k(x−x)  (5)where {circumflex over (x)} is the new (processed) value of a pixel to be computed from the old (unprocessed) value (x), and the local mean (x) of the old value's surrounding and including that pixel. The parameter k is controlled by the ratio of the local variance to the local mean. Dutt V., Greenleaf J. F, “Adaptive speckle reduction filter for log compressed B-scan images”, IEEE Trans. on Medical Imaging, Vol. 15, No. 6, pp. 802-813, December (1996) discloses the same technique and used the same equation in their literature. But they considered statistics of speckles for log compressed ultrasound image the parameter k was chosen as,k=1−{circumflex over (f)}(α)  (6)where {circumflex over (f)} is the statistics and given by,
                              f          ^                =                                            π              2                        24                    ⁢                                                                      D                  ^                                2                            V                        .                                              (        7        )            Here {circumflex over (D)} is an estimate of log compression parameter from the dynamic range and V is the local sample variance.
Jong-Sen Lee, “Speckle Analysis and Smoothing of Synthetic Aperture Radar Images”, Computer Graphics And Image Processing 17, pp. 24-32, (1981) proposed a smoothing algorithm based on local statistics on a fixed window size and was successfully applied to remove speckles form SAR images. They considered multiplicative noise model for speckle where the noise is independent to the signal having mean 1 and variance σv2. The basis of this filter is: in homogeneous region the filtered output is linear average of pixels in the neighborhood, whereas in the region of extremely large intensity variation the output becomes the value of the input pixel itself. The output of the Lee filter is given as,{circumflex over (x)}=x+k(z−v.x).  (8)
Here z is the observed pixel, v=1, and the value of k is calculated as follows:
                              k          =                                    Var              ⁡                              (                x                )                                                                                                          x                    _                                    2                                ⁢                                  σ                  v                  2                                            +                              Var                ⁡                                  (                  x                  )                                                                    ⁢                                  ⁢                              x            _                    =                      z            _                          ⁢                                  ⁢        and        ⁢                                  ⁢                              Var            ⁡                          (              x              )                                =                                                                      Var                  ⁡                                      (                    z                    )                                                  +                                                      z                    _                                    2                                                                              σ                  v                  2                                +                                                      v                    _                                    2                                                      -                                          z                _                            2                                                          (        9        )            
The quantities z and Var(z) are approximated by local mean and local variance of speckle corrupted image.
The main limitation of the Bamber, Dutt, and Lee filters is that the use of too large window introduces a loss of fine details in the image. On the other hand, the use of small window implies insufficient speckle suppression homogeneous region. To avoid this problem adaptive windowing and modified adaptive filtering with variable window size are also proposed in Park J. M., Song W. J., Pearlman W. A., “ Speckle Reduction for SAR Images based on adaptive windowing”, IEE Proceedings Vol. 146, No. 4, August (1999).
Kuan D. T., Sawchuk Alexander A. et al., “Adaptive restoration of images with speckle,” IEEE Trans. Acoustics, Speech and Sig. Proc., Vol. 35, pp. 373-383, March (1987) used same formulation with different assumption of signal model. They assumed that the speckle samples are independent of each other. They derived a local linear minimum mean square (LLMMSE) filter using non-stationary mean and non-stationary variance (NMNV) image model. The correlation properties are also taken into account in their derivation. The parameter k for Kuan filter is determined as,
                    k        =                                            Var              ⁡                              (                x                )                                                                    Var                ⁡                                  (                  x                  )                                            +                              x                _                            +                              Var                ⁡                                  (                  x                  )                                                              .                                    (        10        )            
The MMSE filter proposed by Frost V. S., Stiles J. A., Shanmugan K. S., Holtzman J. C., “A model for radar images and its application to adaptive digital filtering for multiplicative noise,” IEEE Trans. on Pattern Analysis And Machine Intelligence, Vol-4, pp. 157-166, March (1982) is a balance between averaging and all pass filter. The one dimensional impulse response of the MMSE Frost filter is derived as,h(t)=Aαe−α|t|  (11)where A is the normalizing constant and α is the ratio of square root of local variance to local mean of the observed image in a window.
Directional median filter as disclosed in Czerwinski, R. N., Jones, D. L., William D. O'Brien, Jr., “Ultrasound Speckle Reduction by Directional Median Filtering”, Proceedings, International Conference on Image Processing, Vol: 1, (1995) and adaptive weighted median filter as disclosed in Loupas T., McDicken W. N., Allan P. L., “An Adaptive Weighted Median Filter for Speckle Suppression in Medical Ultrasonic Images”, IEEE Trans. on Circuits and Systems, Vol. 36, No. 1, pp. 129-135, January (1989) are also in use for reducing of speckle due to their robustness and edge preserving capability. These filters are nonlinear filters and produce relatively less blurred image. However, their computational complexity is large.
In many cases Maximum-a-posteriori (MAP) filters are used for speckle reduction. MAP filters require assumption about the distribution of the true process and the degradation model. Different MAP estimators are proposed with different assumptions and different complexities as disclosed in Kalaivani S., Narayanan, Wahidabanu R. S. D., “A View on Despeckling in Ultrasound Imaging”, International Journal of Signal Processing, Image Processing and Pattern Recognition Vol. 2, No. 3, pp. 85-97, September (2009).
In Diffusion filtering the nonlinear partial differential equation based smoothing technique utilizing the concept diffusion is proposed by Perona P and Malik J, “Scale-Space and Edge Detection Using Anisotropic Diffusion”, IEEE Trans. on Pattern Analysis And Machine Intelligence, Vol. 4, No.-7, pp. 629-639, July (1990). The diffusion is described by,
                                                        ∂              I                                      ∂              t                                =                      div            ⁡                          [                                                c                  ⁡                                      (                                                                                        ∇                        I                                                                                    )                                                  ⁢                                  ∇                  I                                            ]                                      ⁢                                  ⁢                              I            ⁡                          (                              t                =                0                            )                                =                      I            0                                              (        12        )            where div is the divergence operator and | | is the magnitude, c is the diffusion constant and I0 is the initial image. Two diffusion constants are considered as,
            c      ⁡              (        x        )              =                  1                  1          +                                    (                              x                k                            )                        2                              ⁢                          ⁢      and                  c      ⁡              (        x        )              =          exp      ⁡              (                  -                                    (                              x                k                            )                        2                          )            
In the anisotropic diffusion method, the gradient magnitude is used to detect an image edge or boundary as a step discontinuity in intensity.                If |∇I|>>k then c|∇I|→0, and we have all pass filter,        If |∇I|<<k then c|∇I|→1, and we achieve anisotropic diffusion (Gaussian filtering).        
An edge sensitive diffusion method called speckle reducing anisotropic diffusion (SRAD) has been proposed to suppress speckle while preserving edge information disclosed in Yongjian Yu and Scott T. Acton, “Speckle Reducing Anisotropic Diffusion”, IEEE Trans. on Image Processing, Vol.-11, No.-11, pp. 1260-1270, November (2002). These methods have one common limitation in retaining subtle features such as small cysts and lesions in ultrasound images. A modified SRAD filter, which rely on the Kuan filter rather the Lee filter was developed in Aja-fernandaz S., Alberola-Lopez C., “On the estimation of coefficient of variation for anisotropic diffusion speckle filtering”, IEEE Trans. on Image processing, Vol. 15, No. 9, pp. 2694-2701, September (2005) and this approach is called Detail preserving Anisotropic Diffusion (DPAD). This method is combined with matrix anisotropic diffusion method designed to preserve and enhance small vessel structures referred as oriented speckle reducing anisotropic diffusion disclosed in Krissian K. Fedrij C, “Oriented Speckle reducing anosotropicn diffusion”, IEEE Trans. on Image Processing, Vol. 15, No. 5, pp. 1412-1424, May (2007).
Multiscale methods include wavelet and pyramid based denoising and discussed in several literatures i.e. David L. Donoho, “De-Noising by Soft-Thresholding”, IEEE Trans. on Information Theory, Vol. 41, No. 3, pp. 613-627, May (1995), S. Grace Chang, Bin Yu, Martin Vetterli, “Adaptive Wavelet Thresholding for Image Denoising and Compression” IEEE Trans. on Image Processing, Vol.-9, No.-9, pp. 1532-1546, September (2000), K. P. Soman and K. I. Ramachandran, “Insight into wavelets: From Theory to Practice” PHI (EEE) 2nd Edition, (2005) and Du{hacek over ( )}san Gleich, Mihai Datcu, “Wavelet-Based SAR Image Despeckling and Information Extraction, Using Particle Filter”, IEEE Trans. on Image Processing, Vol. 18, No. 10, pp. 2167-2184, October (2009).
Wavelet denoising attempts to remove whatever noise present and retain whatever signal is present regardless of the frequency content of the signal as mentioned in K. P. Soman and K. I. Ramachandran, “Insight into wavelets: From Theory to Practice” PHI (EEE) 2nd Edition, (2005). It is nothing but shrinkage of wavelet coefficients in wavelet transform domain. Three basic steps are required for wavelet denoising. The steps are as follows:
1. A linear forward wavelet transform,
2. A non-linear shrinking denoising,
3. A linear inverse wavelet transform.
Step1: Calculate the wavelet coefficients of the observed data by applying wavelet transform.
Step2: Modify the detail coefficients to obtain the estimate of the original signal.
Step3: Take the inverse transform of the modified detail coefficient to obtain the denoised signal.
The main challenge of wavelet denoising is the proper choice of shrinkage function and the threshold.
Two categories of thresholding are in use:
Global thresholds: Single threshold (λ) (is chosen to apply globally to all wavelet coefficient Level dependent threshold: Possibly different thresholds are chosen for each resolution level. One should estimate the noise level (σ) to determine the threshold. The above two categories of thresholding include hard thresholding and soft thresholding techniques. The thresholding is discussed in brief, Let w be the observed noisy data, σ the estimated noise level, λ the threshold and Dλ(.) denotes the shrinkage function, which determines how threshold is applied to data. Then modified wavelet coefficients can be given as,
                              w          ^                =                  σ          ·                                    D              λ                        ⁡                          (                              w                σ                            )                                                          (        13        )            
Denoising methods differ in the choices for Dλ(.), λ and σ. Different denoisers consider different shrinkage functions that determine how the threshold is applied, different noise estimates and different shrinkage rules to determine the threshold σ. A few shrinkage functions, which are generally used for denoising, are listed below:
                    Hard        ⁢                                  ⁢        threshold        ⁢                  :                ⁢                                  ⁢                              D            H            λ                    ⁡                      (            w            )                          ⁢                  {                                                                      w                  ,                                                                                                  for                    ⁢                                                                                  ⁢                    all                    ⁢                                                                                  ⁢                                                                w                                                                              >                  λ                                                                                                      0                  ,                                                            otherwise                                                                        (        14        )                                          Soft          ⁢                                          ⁢          threshold          ⁢                      :                    ⁢                                          ⁢                                    D              S              λ                        ⁡                          (              w              )                                      =                              sign            ⁡                          (              w              )                                ⁢                      max            ⁡                          (                              0                ;                                                                          w                                                        -                  λ                                            )                                                          (        15        )                                Garrot        ⁢                  :                ⁢                                  ⁢                              D            G            λ                    ⁡                      (            w            )                          ⁢                  {                                                                                          (                                          w                      -                                                                        λ                          2                                                w                                                              )                                    ,                                                                                                  for                    ⁢                                                                                  ⁢                    all                    ⁢                                                                                  ⁢                                                                w                                                                              >                  λ                                                                                                      0                  ,                                                            otherwise                                                                        (        16        )                                Semisoft        ⁢                  :                ⁢                                  ⁢                              D            SS            λ                    ⁡                      (            w            )                          ⁢                  {                                                                      0                  ,                                                                                                  for                    ⁢                                                                                  ⁢                                                                w                                                                              ≤                                      λ                    1                                                                                                                                                                  sign                      ⁡                                              (                        w                        )                                                              ⁢                                                                                            λ                          2                                                ⁡                                                  (                                                                                                                  w                                                                                      -                                                          λ                              1                                                                                )                                                                                                                      λ                          2                                                -                                                  λ                          1                                                                                                      ,                                                                                                  for                    ⁢                                                                                  ⁢                                          λ                      1                                                        <                                                          w                                                        <                                      λ                    2                                                                                                                        w                  ,                                                                                                  for                    ⁢                                                                                  ⁢                                                                w                                                                              >                                      λ                    2                                                                                                          (        17        )            
The VisuShrink was proposed as a global rule for one-dimensional signals as disclosed in David L. Donoho, “De-Noising by Soft-Thresholding”, IEEE Trans. on Information Theory, Vol. 41, No. 3, pp. 613-627, May (1995). Regardless of the shrinkage function, for a signal size n, with noise from a standard normal distribution N(0,1), the threshold is,λU=√{square root over (2 log(n))}  (18)
If data is not normalized w.r.t noise-standard deviation, first the σ using the equation below is estimated
                              σ          =                                    median              ⁢                              {                                  (                                                                                                                                                                    w                            k                                                                                                    :                        k                                            =                      1                                        ,                    2                    ,                                          3                      ⁢                                                                                          ⁢                      …                      ⁢                                                                                          ⁢                                              n                        2                                                                              )                                }                                      0.6745                          ,                            (        19        )            
VisuShrink is found to yield an overly smoothed estimate. This is because the universal threshold (UT) is derived under the constraint that with high probability the estimate should be at least as smooth as the signal. So the UT tends to be high for large values of n, killing many signal coefficients along with the noise. Thus, the threshold does not adapt well to discontinuities in the signal.
The SureShrink for 1-D data, thresholds derived by minimizing Stein's Unbiased Risk Estimate (SURE) depends on the multiresolution level. It can be generalized to 2-D images in either level or subband dependent manner. The threshold on subband s is:λS=argλ≥0min[SURE(λ, ws)]  (21)
If wavelet decomposition is sparse a hybrid method combining the universal and SURE threshold is preferable over SURE as disclosed in David L. Donoho, “De-Noising by Soft-Thresholding”, IEEE Trans. on Information Theory, Vol. 41, No. 3, pp. 613-627, May (1995). The hybrid method when combined with the soft shrinkage function is referred to as SureShrink in the literature. SureShrink is subband adaptive and a separate threshold is computed for each detail subband. The threshold is also determined for each sub-band assuming a Generalized Gaussian Distribution (GGD) for the wavelet coefficients in each detail sub-band as disclosed in S. Grace Chang, Bin Yu, Martin Vetterli, “Adaptive Wavelet Thresholding for Image Denoising and Compression” IEEE Trans. on Image Processing, Vol.-9, No.-9, pp. 1532-1546, September (2000). The threshold is determined by minimizing Bayes Risk. The pdf associated with GGD is given by,
                                          p            ⁡                          (              x              )                                =                                    [                                                v                  ⁢                                                                          ⁢                                      η                    ⁡                                          (                                              v                        ,                        σ                                            )                                                                                        2                  ⁢                                      Γ                    ⁡                                          (                                              1                        v                                            )                                                                                  ]                        ⁢            exp            ⁢                          {                              -                                  [                                                            η                      ⁡                                              (                                                  v                          ,                          σ                                                )                                                              ⁢                                                                                          x                                                                    v                                                        ]                                            }                                      ⁢                                  ⁢        where        ⁢                                  ⁢                              η            ⁡                          (                              v                ,                σ                            )                                =                                                                      σ                                      -                    1                                                  [                                                      Γ                    ⁡                                          (                                              3                        v                                            )                                                                            Γ                    ⁡                                          (                                              1                        v                                            )                                                                      ]                                            1                2                                      .                                              (        22        )            
Here ν is the shape parameter v describing the exponential rate of decay and σ is the standard deviation. Assuming such a distribution for the wavelet coefficients, we empirically estimate ν and σ for each subband and try to find the threshold T, which minimizes the Bayes Risk as follows:R(T)=E({circumflex over (X)}−X)2=EXEY/X({circumflex over (X)}−X)2  (23)where {circumflex over (X)}=ηT(Y), Y/X˜N(x,σ2) and X˜p(x), ηT(x)=sign(x) max(|x|−T,0)
The optimal threshold T* is given by,
                              T          *                =                  arg          ⁢                                          ⁢                                    min              T                        ⁢                          R              ⁡                              (                T                )                                                                        (        24        )            
T* can be evaluated as,
                              T          *                =                                            σ              ^                        2                                              σ              ^                        X            2                                              (        25        )            where
            σ      ^        =                  median        ⁡                  (                                                Y              ij                                            )                    0.6745        ,Yij ϵSubband HH1, {circumflex over (σ)}X2=√{square root over (max(σY2−{circumflex over (σ)}2,0))} and
      σ    Y    2    =            1              n        2              ⁢                  ∑                  i          ,          j                n            ⁢              Y        ij        2            
BayesShrink performs better than SureShrink in terms of MSE. The reconstruction using BayesShrink is smoother and more visually appealing than the one obtained using SureShrink. Many other types of wavelet thresholding are there in different literature. Two major limitations of transform domain shrinkage methods are that they exhibit (a) pseudo-Gibbs and (b) fake feature types of artifacts in images corrupted with medium to high levels of noise.
Adaptive weighted median filter (AWM): It is already discussed in several literatures that median filters perform better than the linear spatial filters for speckle reduction in ultrasound images. The adaptive weighted median filter proposed by T. Loupas et al is better than the median filter for its edge preserving capability. Their method of calculating the weights of the median filter is based on the manipulating of the local statistics of the image. According to them, the weight co-efficients of the median filter are adjusted by,
                              w          ⁡                      (                          i              ,              j                        )                          =                  [                                    w              ⁡                              (                                                      K                    +                    1                                    ,                                      K                    +                    1                                                  )                                      -                          c              ⁢                                                          ⁢              d              ⁢                                                σ                  2                                m                                              ]                                    (        28        )            where c is the scaling constant, m, σ2 are the local mean and variance inside the 2K +1 by 2 K +1 window, d is the distance of the point (i, j) from the center of the window (K +1), K +1) and [x] returns the nearest integer to x.
The problem of this method is the selection of the constant c and the center weight w(K+1, K+1), which influence the results. No definite method is given to determine the value of these quantities.
US20070071292 discloses image processing adapts to speckle. Speckle is identified from signal transitions. For example, peaks, valleys or mean crossings of image signals as a function of space or spatial location are identified. A speckle characteristic, such as speckle size, is estimated from the signal transitions. The estimation may be limited to soft tissue regions to reduce the effects of specular targets and noise on speckle estimation. The speckle is estimated for local regions or an entire image. By estimating speckle for local regions, image processing may account adaptively for regional variation in speckle size.
US20090240144 discloses systems and methods for suppressing speckle noise in ultrasound imaging. In an embodiment, speckle noise suppression is provided by incoherently summing echo waves that impinge the active aperture of the transducers. This incoherent summation prevents echo waves from destructively interfering and therefore prevents the signal ‘nulls’ that characterize speckle noise. In an exemplary embodiment, the incoherent summation is performed by sub-dividing a transducer into a plurality of smaller transducers and incoherently summing the electrical signals from the smaller transducers. In one exemplary embodiment, each of the smaller transducers is coupled to a separate rectifier, which rectifies the electrical signal from the respective transducer into a rectified signal. The rectified signals from the rectifiers are then summed to provide the incoherent summation.
U.S. Pat. No. 5,653,235 discloses a system and a method for generating an ultrasound image of an interrogation region in an object with a transducer with a two-dimensional array of transducer elements, includes the steps of generating an ultrasound beam by activating many transducer elements of the two-dimensional array; electronically controlling the beam to illuminate substantially the same region from at least two orientations; capturing the echoes generated with the beam illuminating the object at different orientations; and analyzing the echoes from all directions to produce an image of the region of the object. The aperture of the transducer generating the beam is at least substantially equal to the aperture generated by a linear array of transducer elements extending across the substantially shortest distance between two opposite edges on the two-dimensional array. The ultrasound power emitted from the transducer elements is not spatially uniform, and the multiple echoes reduce speckle in the image.
U.S. Pat. No. 6,517,486 discloses a compounding method for reducing speckle noise applied in an ultrasound imaging apparatus is disclosed. The compounding method includes the steps of providing an object, measuring the object for obtaining a reference image by the ultrasound imaging apparatus, applying an external force to the object to deform the object, measuring the deformed object for obtaining an deformed object image at the same position, estimating an in-plane displacement field of the deformed object image for correcting an in-plane motion of the object to obtain a corrected image, and compounding the reference image with the corrected image to obtain a compounded image of the object for achieving the speckle noise reduction.
U.S. Pat. No. 5,497,777 discloses the enhancement of ultrasound images is provided through the filtering of signal dependent noise such as speckle noise by dividing the signal into selective subintervals and utilizing discrete wavelet transform and the identification and selection of those wavelet transform coefficients primarily including signal and not those primarily including signal dependent noise.
U.S. Pat. No. 5,409,007 discloses a method for reducing speckle artifact in an ultrasound image using a two-dimensional median filter having a diamond-shaped five-point kernel. The entire pixel image data is passed through the filter in a manner such that the center point of the kernel is effectively stepped down each range vector in sequence. The magnitudes of the pixel data at each of the five points in the kernel are compared and the value which has the middle magnitude is adopted as a new pixel value, which is substituted for the old pixel value at the center point. After a new filtered vector has been formed from the new pixel values produced at successive center points by stepping down one acoustic vector, the kernel is shifted by one vector and stepped down range again. This process continues through the entire set of vectors until a new set of filtered vectors is formed. This filter will remove speckle holes on the order of one pixel in size while preserving good edge definition.
U.S. Pat. No. 6,454,715 discloses a methods and apparatus for blood speckle detection for enhanced intravascular ultrasound imaging. The present invention utilizes the fact that the energy scattering strength from blood exhibits a high frequency dependency, while the scattering strength from tissue lacks a strong frequency dependency. In specific embodiments, the present invention may provide a particularly simple and useful solution for addressing the problem of blood speckle in intravascular ultrasound imaging, especially in situations where the blood may have a scattering strength similar to that of tissue and/or where the blood is moving slowly or not at all.
US 20040127795 discloses a method and apparatus for smoothing speckle pattern and increasing contrast resolution in ultrasound images is provided. Compared to other frequency compounding techniques, wide-band harmonic frequency compounding reduces speckle noise without sacrificing the resolution. Compared to spatial compounding, wide-band harmonic frequency compounding is more robust against tissue motion because sequential vectors rather than frames are summed together for compounding. The method and apparatus is implemented by transmitting two or more firings, combining two or more of the firings coherently to extract the tissue-generated harmonic components, detecting the outputs of the coherent sums and detecting one or more firings before coherent sum, and finally combining all detected outputs to form the compounded image. The method and apparatus sums wide-band fundamental and wide-band harmonic images after detection to form a compounded image. Unlike other frequency compounding methods, both transmit and receive signals are wide-band and no narrow-band filters are necessary. Multiple firings with two or more different transmit waveforms are transmitted to each focal zone.
US 20070065009 discloses a method for enhancing an ultrasound image is provided, wherein the ultrasound image is segmented into a feature region and a non-feature region, while sufficiently utilizing features contained in the ultrasound image, in particular including some inconspicuous features. The enhanced image according to present invention is not susceptive of the image segmentation and avoid dependence of the enhancement effect on the segmentation template, so as not to produce an evident artificial boundary between the feature region and the non-feature region but to highlight some special information in the image and to remove or mitigate invalid information. Thus the enhanced ultrasound image is particularly suitable for the visual system of the human beings.
US 20080181476 discloses a methods and systems for enhancing an image exhibiting speckle noise are provided. An image exhibiting the speckle noise is received and a coefficient of variation is estimated in a part of the received image. Either a detail tuning parameter or a smooth tuning parameter are selected based on the estimated coefficient of variation. A maximum likelihood (ML) filter is configured with the selected tuning parameter and the configured ML filter is applied to the part of the received image.
U.S. Pat. No. 7,720,266 discloses a method for enhancing an ultrasound image is provided, wherein the ultrasound image is segmented into a feature region and a non-feature region, while sufficiently utilizing features contained in the ultrasound image, in particular including some inconspicuous features. The enhanced image according to present invention is not susceptive of the image segmentation and avoid dependence of the enhancement effect on the segmentation template, so as not to produce an evident artificial boundary between the feature region and the non-feature region but to highlight some special information in the image and to remove or mitigate invalid information. Thus the enhanced ultrasound image is particularly suitable for the visual system of the human beings.
The drawbacks of the above mentioned prior art is that the amount of data to be handled during speckle reduction is large and thus increases complexity. Furthermore, after noise reduction from raw data, interpolation is to be done for scan conversion, which further increases complexity. The Conventional interpolation is a complex technique for the reduction of speckle in ultra sound imaging. The computational complexity and the number of building blocks are more in case of the above mentioned prior art. The interpolation stage in the prior art increases the loss of information and does not provide better output image quality. The prior art systems and method cannot be used with SR (super resolution) ultrasound image reconstruction techniques.
Thus there is a need to provide an improved ultrasound imaging method/technique for speckle reduction/suppression in an ultra sound imaging system that does not require any conventional interpolation during scan conversion and can be also implemented with Ultrasound SR reconstruction technique from polar format data. Further to provide better image quality by removing speckles and preserving edges.
Further reducing the overall computational complexity and the number of building blocks of conventional ultrasound imaging system and it can be used with SR ultrasound image reconstruction techniques.
Scan conversion means evaluating the pixel values at the grid points in rectangular co-ordinate system. In ultrasound sector scanner the pixel values are available at the grid points in polar co-ordinate system after scanning the object. Conventionally, in scan conversion process the pixel values at the grid points in rectangular co-ordinate system are evaluated by using interpolation techniques using the pixel values available in polar co-ordinate system. This scan conversion makes it possible to display the ultrasound image in video monitor which supports inputs in rectangular co-ordinate system only to display the image. After scan conversion the speckle reduction technique is generally employed to remove the speckle noise from the scan converted image. In the present technique we are avoiding this interpolation.
In conventional techniques, the pixel values at Q are calculated by using a suitable interpolation technique. Here, the values of P are first calculated along radial direction using interpolation technique. After calculation of the pixel values at P, the pixel values at Q are calculated using interpolation technique. This completes the scan conversion process. Speckle reduction is applied after this scan conversion process i.e. the pixel values at the grid points (Q points) in the rectangular co-ordinate points are re-calculated by using the speckle reduction filtering algorithms (e.g. Lee, Kuan, median, weighted median, adaptive weighted median filters etc.). The drawback of this conventional procedure is that the interpolation in the scan conversion stage makes the noise more coloured and the effect of noise is spread from a smaller region to relatively bigger region. Moreover, we lose some information in the scan conversion process due to low pass nature of interpolation operation. This degrades the performance of the conventional procedure. In the present technique, we reduce the loss of information in the in scan conversion process since scan conversion is performed through filtering.
The inventors have found an improved speckle reduction method where the prior speckle reduction techniques can be used to obtain better quality of output image. The inventors have used noise/speckle reduction filter during scan conversion instead of using them after scan conversion. Further the inventors also proposed an improved method for speckle reduction using an improved filter which gives better quality of image if it is applied in old speckle reduction technique. However applying it the improved speckle reduction technique/method provides much better quality of image than that of old speckle reduction technique.